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Abstract 

This paper deals with a hyperbolie system of two nonlinear eonservation laws, 
where the phase spaee eontains two eontaet manifolds. The governing equations 
are modelling bidisperse suspensions, whieh eonsist of two types of small parti- 
eles that are dispersed in a viseous fluid and differ in size and viseosity. For eertain 
parameter ehoiees quasi-umbilie points and a eontaet manifold in the interior of 
the phase spaee are deteeted. The dependanee of the solutions strueture on this 
eontaet manifold is examined. The elementary waves that start in the origin of the 
phase space are classified. Prototypic Riemann problems that connect the origin 
to any point in the state space and that connect any state in the state space to the 
maximum line are solved semi-analytically. 
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1. Introduction 

Polydisperse suspensions can be described by balance equations as N super¬ 
imposed continuous phases, where particles of species i, associated with a volume 
fraction 0*, distinguish in properties like size, density and viscosity [|7l[I3; mod¬ 
els with similar solution structure describe traffic and pedestrian flows [|4l[8l. For 
the considered model of bidisperse suspension, the solution structure of the solu¬ 
tion to the initial value problem for standard batch settling tests has been studied 
for the cases when strict hyperbolicity is assured [|5l| and when the phase space 
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provides elliptic regions B . The focus of this contribution is on the impact of a 
contact manifold in the interior of the phase space, that emerges for certain pa¬ 
rameter settings. This contact manifold has the physical property of coinciding 
particle settling velocities, and thus has a practical relevance, since one general 
goal in the process control of solid-liquid separation processes is to reduce segre¬ 
gation effects [fTOl . 

The generic form of kinematic sedimentation models for poly disperse suspen¬ 
sions consists of the system of N first-order hyperbolic equations 

dt<Pi + dJ\{<^) = 0, /,($) = t = ( 1 ) 


where t is time, x is depth and the velocity components Ui(<h) depend on the 
concentration vector $ = (^i,..., The unknown $ denotes the vector 

of volume fractions of the solids phases and is contained by the phase space of 
physically relevant concentrations 


V 


$co 


— ( 01 , . . . , £ 


pAf . 


01 > 0, , 0Ar >0, 1 

0 := 01 -f . . . -f 0Ar < 0“ J ’ 


( 2 ) 


where the total concentration 0 := 0i -|- ... -f 0Ar is bounded from above by the 
maximum packing concentration 0“. The maximum packing manifold 


d°° {$ — (01,..., 0iv)'^ : 0 01 -f • • • -f 0Ar — 0°°}. (3) 


is the set of all maximal states. 

The resulting system of conservation laws is actually a system of mass bal¬ 
ances for different solids species, where the nonlinear flux function 

/(«):= (/i(1>),...,/iv(1>))’' 

can be derived from the corresponding momentum balances [|71[I^|24l| The com¬ 
ponents describe the flow process of the dispersed solids phases in a liquid, where 
the dispersed phases are considered as a continuum. The flux function has com¬ 
ponents 

/i($) = 0iUi(<h), Ui(<h) = Mi(<h) - i = (4) 

with u = (mi($), ... ,M 7 v(*h))^, where the absolute velocity Vi = Ui($) of a 
representative solids particle depends on a linear combination of the solid-fluid 
relative (“slip”) velocities 

Mi($) := t;i($)-Uf, (5) 
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which are relative to the fluid veloeity Uf. The flux funetion model Q is elosed 
by speeifying the relative veloeity as 

Mi($) = VooiVi{^), (6) 


where the constant Vad is the Stokes veloeity, whieh quantifies the settling veloeity 
of a single partiele in a fluid, and 14(<h) is the hindered-settling veloeity that is an 
non-inereasing funetion of the eomponents of <h, see O. Following Riehardson 
and Zaki ||2^ . the hindered-settling veloeity is set as 





if 0 < 0 < 
otherwise. 


i = 1, 


.N, 


(7) 


where the exponent n* > 1 aeeounts for the slow down of the proeess at inereasing 
eoneentrations. Assumptions @ and Q ean be eombined as 

= t^oo*(l - 0)”'"' (8) 


for <h G . 

Strietly hyperbolie systems of eonservation laws, where the eigenvalues of 
the Jaeobian matrix of the flux funetion are real and distinet, provide a relatively 
well understood framework for the solution of Riemann problems ffTSl . In [[T2|, 
striet hyperbolieity of the system ([T]) with flux function Q has been first shown 
for N = 2 and later on in 0 for general N, but up to then only for coineiding 
hindered-settling faetors 14(0) = 14(0) = ■ ■ ■ = 14v(0), whieh depend on the 
total eoneentration 0. For a model with the more general hindranee faetor (|7]), 
this implies to have eonstant exponents ni = ^2 = ■ ■ ■ = tiat, a restrietion that 
turned out to be unneeessary: In [[3l, it is shown that striet hyperbolieity also holds 
for general N > 2 and relative veloeities of form 14($) = 14(0) as long as the 
inequality 

M'(l-0)-Mi<O (9) 

holds for alH = 1,..., A^, where u' denotes the derivative with respeet to 0; 
in this situation the only restrietion on the hindered-settling funetion 14(*h) = 
14(0) is to depend on the total eoneentration 0. The inequality @ is satisfied in 
partieular when the relative veloeities are ordered as ui > U 2 > ■ ■ ■ > «„ for any 
0. This holds in the ease of the hindered-settling funetion 0 for the situation if 
the parameters are ordered as 


^^ooi > > ■■■ > VocN with ni <n 2 < ■■■ < n^- (10) 
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In [fTTlI . a secular equation framework was established that allows to verify strict 
hyperbolicity by checking a simple algebraic criterion; this framework has been 
applied to the considered model with a Richardson-Zaki hindered settling func¬ 
tion having constant exponent. In this framework was adapted to the same 
general model setting as in [|3l, i.e. with size-dependent hindered settling factors 
that not necessarily take the form Q. Subsequently, in IfTTl the secular equation 
framework was applied to a series of choices of hindered-settling functions. 

Preliminary numerical simulations of Riemann problems for N = 3 with arbi¬ 
trary parameter choices (not exposed here) showed that strict hyperbolicity might 
fail as coincidence of eigenvalues occurs, providing an abrupt change of the so¬ 
lution structure. The fact that this phenomenon of abrupt change already appears 
for N = 2 made us to look for analytical insights in this situation, which lead to 
the present contribution. 

This contribution deals with the wave classification for 2 x 2 systems of conser¬ 
vation laws that arise as one-dimensional kinematic models for the sedimentation 
of bidisperse suspensions. The analytical examination of bidisperse suspensions 
gives insights to flow properties of polydisperse suspensions, which are mixtures 
of small solid particles dispersed in a viscous fluid. In this contribution the focus 
is on models for particle suspension where all particles are assumed to have the 
same density. Specifically, in this contribution, the properties of the 2x2 system 
([T]) (with N = 2), flux function Q and closures Q are studied. The model of 
our interest contemplates the following specifications, which is done in opposition 
to ( fTO] ), which would guarantee strict hyperbolicity: 


(51) Wool > Voo2 > 0, 

(52) rii > n2 > 1, 

(53) 0“ = 1. 

For N = 2 and Vi{(j)) given by Q the flux function /(<T) = (/i(*I’), 
takes the form 


/i($) = 0i(^nooi(l - - 0)”^ ^ - noo202(l - (/>)”" 

/2($) = (j)2(v^2{l - 02)(1 - 0 )'""“^ - nool 0 l(l “ , 


for values <I) G V^oo and/i($) = f 2 {^) = 0 otherwise. A special interest consists 
in the classification of the solution structure of the Riemann problem 



X < 0, 
X > 0. 


(II) 
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For convenience, the Riemann problem consisting of the system of PDFs ([T]) with 
initial condition ( [TT] ) is referred to as RP(<h“, $+), with left and right values 
and respectively, to be specified. 

The application of this model is the batch settling process of an initially homo¬ 
geneous suspension in a closed container described by the initial-boundary value 
problem 


dt(f>i + = 0, i = l,2, 

<h(0, x) = <ho(3^)5 0 < X < L, (12) 

/,($)= 0, xe{0,L}, * = 1,2, (13) 


where L is the domain height and the components of the 
fi^) = are given by Because of the 

condition ( fT3| ), the initial-boundary data ( [T^ and ( [T3] ) can 
Cauchy data 


flux-density vector 
zero-flux boundary 
be replaced by the 


<l)(0,x) = <ho(a;) 


{ O for X < 0, 

$0 for 0 < X < L, 
for X > L, 


(14) 


where O := (0, 0)^ is the origin and is a state on the maximum concentration 
manifold Q. Therefore, the Riemann problems RP(0, <!)) and RP($, <1)°°) are of 
particular interest. This contribution reveals analytical insights into the solution 
structure of the Riemann problem RP(0, <!)). From an application point of view, 
the Riemann problem RP(0, <!)) describes the interactions on the upper interface 
between clear liquid and initially homogeneous suspension during a batch settling 
process. 

With respect to related work, several studies on weakly hyperbolic systems, 
i.e. systems that are hyperbolic but not strictly hyperbolic, are developed for mod¬ 
els of multi-phase flow in porous media. For three-phase flow in porous media, 
the Corey model with convex permeability leads to a single isolated point, the 
so-called umbilic point, in which strict hyperbolicity fails [[I3l [20l |21 HSl l29ll : 
another loss of hyperbolicity occurs when a the phase space contains an elliptic 
region as it occurs for the Stone model [flSll . To solve Riemann problems, the 
wave-curve method has been applied, in which a sequence of elementary waves 
are connected. When strict hyperbolicity fails, it is not sufficient to consider the 
method of Liu to deal with non-convex fluxes; rather one has also to considers 
non-local branches of the Hugoniot locus [lT4l . The wave-curve method has been 
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applied to the injection problem, where a gas-water mixture is injected in a porous 
medium containing oil [|2|. For a system of two conservation laws with a quadratic 
flux function the solution in the neighborhood of the umbilic point has been clas¬ 
sified in ifT^ 1^ 1^ . Following the idea of studying the solution behavior in the 
neighborhood of an umbilic point, in the case of a quadratic flux function four 
different types of umbilic points corresponding to different shapes of the close by 
integral curves could be identified [|29l . In the Corey model with convex perme¬ 
ability only two types of umbilic points occur [l25ll . According to the proposed 
classification, certain types of Riemann problems have been considered in fl20l . 
A systematic classification of solutions of the Riemann problem for non-strictly 
hyperbolic systems of two conservations laws, which count with an umbilic point 
and with the identity viscosity matrix has been carried out in [|^ [311 . A non¬ 
local Hugoniot locus leads to non-classical waves and in some cases to transitional 
shocks idl |23l ■ These transitional shocks are sensitive to the regularization by a 
non-identical viscosity matrix. 


2. Basic definitions 


In this paragraph several definitions lIT^ [T^ l26ll are collocated in order to 
facilitate the appropriate classification for the system under study. 

Definition 1. The Hugoniot locus of a state denoted as TLif^~), is the set of 
all states $■*“ that satisfy the Rankine-Hugoniot condition 

/(^>+)-/($-)=a($+-$-), (15) 


where a = (t(<I) , <h+) is the propagation velocity of the discontinuity. 

The shock classification according to Lax ll23 is used in order to refer to 
a subset of the Hugoniot locus that corresponds to a certain wave family. The 
corresponding admissible shocks are classified depending on inequalities between 
the first and the second eigenvalues at both sides of the discontinuity and the 
discontinuity speed itself, see e.g. [l22l[30l . 


Definition 2. Three kinds of classical admissible shocks can be distinguished. 
The classification applies between a left state and a right state $''■ which 
are connected by the Rankine-Hugoniot condition ( [T?] ) with jump velocity a = 
<!)+).• 

1- Lax shock: Ai(<I>’'') <a< Ai(<h“) and a < A 2 ($’'') 

2- Lax shock: A 2 (‘h^) < cr < A 2 (*h“) and Ai(<h“) < a 
Over-compressive shock (OC): A 2 ($''') < a < Ai($“) 
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Left- and right-characteristic shocks are included in this shock type definition. 
They occur when a shock speed coincides with the characteristic speed. Whereas 
by definition an over-compressible shock cannot be characteristic, the limit of the 
inequalities above are included in the 1-Lax and 2-Lax shocks, which are also 
called first and second shock waves. 

An inflection manifold Xj is determined for all states where the Lth eigen¬ 
value attains a maximum or minimum value along the integral curve of the same 
family. 

Definition 3 (Inflection curve). The i-th inflection manifold is defined as 

X, := {$ G X>$oo : VAi(<h) ■ ri(<h) = O}, 

where Aj is the i-th eigenvalue of the Jacobian matrix of the flux function and ri is 
the corresponding eigenvector. 

In the sense of the invariant manifolds defined in [[3^ . we introduce the fol¬ 
lowing concept 

Definition 4. An i-th contact manifold occurs when the i-th integral curve passing 
through a state <I>° coincides with a part of the Hugoniot locus such that 

any state on this intersection satisfies 

A,($'’) = a(<h°, <h) = A,($), 

for the shock speed ct(<I)°, <I)). 

A necessary condition for establishing an i-th contact manifold is that the in¬ 
tegral curve is not a rarefaction in the usual sense, but that the characteristic speed 
is fixed along the curve. 

For states on a contact manifold the following transitivity rule holds. If and 
^2 mutually belong to the Hugoniot locus of the other, connected by a shock of 
speed a = a(<hi, < 1 ) 2 ) and if $2 is on the Hugoniot locus of a state $3 by a shock 
of the same speed a, then belongs to and cr(<hi, $3) also coincides 

with a. This is the essence of the Triple Shock Rule [IT3l [T4l [T9l . Another useful 
version establish the following. 

Lemma 1. Let $ 1 , <1)25 *^3 be non-collinear states such that <hi, < 1)2 belong to 
'^(<1)3) and <I)i belongs to TL{fb2), then (t(<I)i, <^2) = cr{^2: ‘^’3) = <^(‘^’15 ‘^*3)- 

On a contact manifold, rarefactions are indistinguishable from shocks; all 
speeds match a characteristic speed. Remarkably, a Hugoniot locus with constant 
characteristic speed is planar and coincides with the integral curves [|3^ . 


7 



3. Contact manifold 

If the specifications (SI) and (S2) of the considered model hold, then a contact 
manifold inside the phase space can be identified. This manifold turns out to be 
decisive for the characterization of solutions of Riemann problems, in particular 
because the origin is connected to this manifold by a right characteristic shock. 

Definition 5. A set C(<I)*) is defined as a subset of the phase space which 
contains a state = ( 0 ^, fif) that is connected to other states by the property 

C($*) := {$ e Dd,- : Ui(<l>) = r;2(<^’) = (16) 

The state is a representative of the set C(<l>*). The definition of a set C(<l>*) 
unifies two complementary properties: 

1 . yi($) = y 2 ($) for all $ G C(<I)*), 

2 . Ui($“) = Ui($+) for <I)+ G C(<I)*) and z = 1, 2. 

Property (1) describes the local coincidence of velocities of different phases in 
one particular state, whereas property (2) describes the constancy of the velocities 
of a single family along the manifold. 

The following Lemma is a generic result, stating that, a set C (<!)*) is a contact 
manifold whenever the flux function has a certain structure. 

Lemma 2. If the flux function of the system has the structure 

/,($) = (IV) 

then, for any G the set C(<I>*) is a contact manifold with constant shock 

speed 

a($-,$+) = r;i(<l>*). (18) 

Moreover, if the structure of the flux function 0 is such that the absolute velocity 
V depends on the total concentration cj), namely 

r;i(<l>) = vffl), (19) 

then the contact manifold C(<1>*) consists of the line 

C($*) = {$ G : 01 + 02 = 0*}, (20) 

where 0 * = 0 ^ + 02 is the total concentration of $*. 
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Proof. By definition ( [T^ any states $ , $+ G C($*) satisfy that 

Wi(<h-) = ^;2(<h-) = ^;i($+) = ^;2($+) = 

sueh that one ean faetorize 


0+^,(<h+) - = ^i(<h*)(0+ - 0-), z = 1, 2 . ( 21 ) 

From the speeifie strueture of the flux funetion ( [T7] ) one reeognizes that equation 
( [2T] ) states the Rankine-Hugoniot eondition ( [T5] ) with shoek speed ( [TS] ). This es¬ 
tablishes that all states on a set C(<1)*) belong mutually to the Hugoniot loeus of 
eaeh other and any two states on a set C(<h*) ean be eonneeted by a shoek of the 
same speed. 

Sinee shoek speeds loeally eonverge to an eigenvalue, on a manifold with 
eonstant shoek speed any two states and have an eigenvalue that eoineides 
with the shoek speed, 

A,(<h-) = A,($+)=ni(<h*), 


sueh that Definition H] of a eontaet manifold follows. 

The shape of a manifold C($*) to be a line ean be dedueed from property 
( fT9l ) whieh assures that the absolute veloeities are eonstant on the line ( [20l ), i.e. 
ui($) = U 2 (<h) for all $ G C(<1)*). □ 


With this Lemma, nothing yet is said about the existenee of sueh a eontaet 
manifold for the eonsidered model equation. Neither it is deeided to whieh ehar- 
aeteristie family the eontaet manifold belongs, whether to the first or the seeond 
one. 

A further eharaeterization of a set C($*) is developed in the sequel, starting 
from properties that ean be derived from the generie strueture of the model, lead¬ 
ing to properties that depend on partieular model speeifieations. A property that 
ean be used in several instanees is 


r;i($) = r;j(<l)) Ui{^) = i,j G {1, 2}, (22) 

whieh follows direetly from @. Property ( |22l ) assures that eondition ( [T9| ) is sat¬ 
isfied for the model under eonsideration, where the flux funetion has the strueture 
Q eomplemented by eonstitutive assumptions @ and (|7]). 

Up to now it has been shown that any set C(<h*) is a eontaet manifold, whieh 
has in addition the shape of a line. In the next Lemma it is shown that for the 
eonsidered model speeifieations sueh manifolds effeetively exist. 
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Lemma 3. If conditions (SI) and (S2) are satisfied, then two distinct contact man¬ 
ifolds exist in the domain P$oo, namely 

C($“) := = = 01 + 02 = 0“}, (23) 

C($") := {$ = (01, 02)^ : 01 + 02 = 0"}. (24) 

The manifold C($°°) is represented by any state G d°°. A representative state 
that defines the contact manifold C(<l>^) that is distinct to the manifold C($°°) 
can be identified as 

= (0^ 0)^, 0" = 1 - (25) 


Proof First, it is shown that the boundary d°° is a contact manifold. For any 
state G d°° one has Mi(<h°°) = M 2 (<h“) = 0. By property ( |22l ) one gets 
fi(<h°°) = r; 2 (‘h°°) = 0 for any state G d°°, satisfying the Definition ( [T^ of 
the set C(<h*), which by Lemmais a contact manifold. 

To show that there is an additional contact manifold C ($""), one has to find a 
set of states ^ C($“) such that fi(<h^) = n 2 (*^*'')- Because of property ( |22l ) it 
is equivalent to find a such that Mi(<h"') = M 2 (<h^). Resolving 


Vloo(l - = t^2oo(l - 0"") 


x\n2 


with respect to (fi^, which is the total concentration of any representant gives 
( |25] ). The conditions on the parameters (SI) and (S2) guarantees that ((A G (0,1) 
exists. □ 


Throughout this work, the notation C(<h*) is used to refer to a generic contact 
manifold, whereas C(<h'') and C(<h°°) refer to specific contact manifolds with as¬ 
signed representative values and <h“, respectively. It turns out that the contact 
manifolds C(<I)'") and C(<I)“) form transversal branches of the Hugoniot locus of 
the origin. 

Lemma 4 (Hugoniot locus of origin). The Hugoniot locus 'H(O) of the origin 
O = (0, 0)^ consists of four branches: the two coordinate axes as local branches, 

:= {<h = ( 0 , 0)T, 0 G [0,1]}, := {<h = (0, 0 )^, 0 G [0,1]}, (26) 


with variable shock speed 

cr(0, $) = nj(<h) for ^ E d\ i = l,2, 


(27) 
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and the two contact manifolds and identified as ( [2^ and ( |24l ), as 

transversal branches with constant shock speed 


a(0,$) = ui($) =V 2 ($) for (28) 

Proof With the structure of the flux function Q the Rankine-Hugoniot condition 
( [T5] ) connecting the origin O with any state $ takes the form 

a{p,^)(j)i = Vi{^)(j)u ? = 1 , 2 . 


This system of two equations has two possible kinds of solution: On the local 
branches on the axes, and 8 “^, one of the two equations becomes obsolete, since 
both sides vanish on the considered axis. Therefore, the velocity of the remaining 
equation determines the shock speed as p7] ). On the transversal branches, C(<f)^) 
and C(<f>°°), no such cancellation occurs such that for a solution the velocities are 
required to be equal, giving speed (|2^. □ 


The triple shock rule as stated in Lemma [^applies to the connection of the origin 
to any state $ G C(<f>*) having speed ( [2^ with any middle state G C'($*) 
having speed ([T^. Indeed, O, <I)^ are not collinear, thus we have ct(0, <I)^) = 
, $) = a{0, $). This means that any (shock) solution O —)■ $ can be 
constructed by the shock O —)■ followed by a second shock —)■ $ of same 
speed; both solutions determines the same wave pattern, so the same solution. 

Another useful property is the convertibility of the relative velocity «,($) with 
the absolute velocity ni($) on the edges: 


Lemma 5. For states on the edges, $ G 0*, i = 1,2, one has 




(29) 


Proof A state $ G 0* on an axis has the representation <I) = + 0 ^^ 262 , 

where e^, fc = 1,2, are unit basic vectors and 5ik is the Kronecker symbol. 

Then, the definition of the relative velocity t'i(<I>) reduces to 

ni(<I)) = Mi($) - = Mi($) - (j)Ui{^) = (1 - (j))Ui{^), 


which establishes the announcement (29). 


□ 
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4. Characteristic speeds 


System ([T]) can be written in quasi-linear form as 

+= 0 , 


where J($) is the Jacobian matrix of the vector valued flux function /($) ;= 
(/i($),/ at($))"'". The structure of the Jacobian matrix is examined for gen¬ 
eral N in [|3l|TT]|; for iV = 2 it becomes 


J($) = 


(Jiim 


J22($V 


fvi+ M1101 
\ U2i4>2 


M1201 
V2 + U22<p2 


(30) 


or, componentwise, 


Jij Vi6ij -f cj)i'Uij^ ijj 1,2, 

where 5ij is the Kronecker symbol, Vi is the absolute velocity @, and Uij is spec¬ 
ified as 


Uij = Uij{^) = m '($) - i, j = 1,2, (31) 

where It'($)= ^^($) U2{^)j . 

Recall that a system is strictly hyperbolic if the Jacobian matrix of the flux 
function has distinct real eigenvalues. For N = 2, a system of conservation laws 
([T]) is strictly hyperbolic if the discriminant 

:= (Jn($) - J22($))' - 4 Ji 2($) J2i(<l>) 

of the Jacobian matrix ( [30l ) of the flux function / is positive. For a hindered 
settling factor given by Q strict hyperbolicity holds for identical exponents rii = 
n 2 , see [fT^ . which is proofed by showing algebraically that A<j, > 0 for <1) in the 
interior of the phase space for a specification (S3) with 0“ = 1; moreover, 
strict hyperbolicity also holds for the case with different exponents ni ^ n 2 - This 
can be shown by a straightforward calculation [l3]|, which yields the discriminant 
composed by a sum of a square and a positive term 

A$ = [(ni0i - 1 )mi(<F) - (n2(/>2 - 1 )m 2(*1))]^ -f 4nin20i02Mi(*l')M2(*l)). (32) 

This term is positive because of conditions (SI) and (S2) together with the bounds 
01 , 02 A 0, 0 = 01 -f 02 < 1 given by definition of the invariance domain Vq,<x, 
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and the definition of «,($). The positivity of the diseriminant A$ > 0 indieates 
striet hyperbolieity of the system ([T]). 

In the sequel, the strueture of the eigensystem of the Jaeobian matrix J of the 
flux funetion Q is derived. The eigenvalue ealeulation is tedious but straightfor¬ 
ward, leading to the following Lemma. 


Lemma 6. The eigenvalues of the Jacobian matrix ( 301 are calculated as 


Ai,2 = ^ ^ [nifiui + n2(t)2U2\ ± ^ 


where A$ is the discriminant that has the form (32). 


With help of this eigenvalue speeifieation, the eharaeterization of the eontaet 
manifold ean be eompleted by the following Theorem, whieh is formulated for 
any set (3(<I)*), sueh that it applies partieularly for the sets and C(<I)“). 


Theorem 1. The eigenvalues for any state <I) G (3(<I)*) can be specified as 


Ai($*) = and A 2 ($) = ni($*). (33) 


Therefore, any setC{^*) is a contact manifold with respect to the second charac¬ 
teristic family. 

Proof Sinee any state $ in the set C(<I)*) is (by definition) eharaeterized by the 
property vi(<h) = = 'yi(‘^**), the eigenvalues beeome 

Ai,2(*^>) = ^ [nifiiufi^) + n24>2U2{^)] ± ^ 

= vi(<l)*) - 2^/^ 2 

giving the pair of eigenvalues ( [3^ . Here, the last step is justified by property ( |22l ), 
i.e. the equivalence of the equalities ni($) = n 2 (‘^’) and Mi(<h) = M 2 ($): Namely, 
one can relate the term in the parenthesis to the discriminant ( |32l ) as 

[riifiiUi -1- 71202^2]^ = [nifilUi - n 2 (() 2 U 2 f + 4nin20102MlM2 = A$, 


which is valid in the special case when mi(<1 >) = M 2 (*h). 

The association of the contact manifold to the second family can be seen by the 
fact that the eigenvalue A 2 , according to the established values in ([3^, is constant 
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on C($*), A 2 ($) = ni($*) for all $ G With the shock speed established 

in ( [TS] ) in Lemmafor the connection of any two states on C($*), one gets 

A2 ($-) =cr($-,<f)+) = A2(<f>+) 


for all $ , <f)+ G C($*) establishing that any set C(<f)*) is a contact manifold. □ 


Theoremapplies to both contact manifold C(<f)^) and For the line 

C(<f)°°) the contact manifold is in addition characteristic with respect to the first 
characteristic family. 

For <f) G C(<F^) the smaller eigenvalue Ai depends on the discriminant and only 
the bigger eigenvalue A 2 (‘F) = is constant, Thus, a contact manifold C($*) 

is an integral curve of the second family. On the maximum packing manifold 
for G 9“ = C(<f)“), where 01 + 02 = 1 holds, both eigenvalues vanish: 


Ai($) = A2(<f>) =0. 

Theorem [^states explicit expressions for the eigenvalues on the contact man¬ 
ifold and which form part of the Hugoniot locus of the origin as 

identified in Lemma The remaining eigenvalues along the Hugoniot locus of 
the origin can be evaluated directly from the general eigenvalues according to 
Lemma For on an edge d\ i = 1, 2, i.e., 0, = 0 for some f = 1, 2, as 
specified in ( | 2 ^ , the discriminant reduces to A$ = [M 3 _i(<f))(n 3 _j 03 _j — 1 )]^, 
where the subindex 3 — z becomes 2 on the axis and 1 on the axis d"^. A simple 
calculation shows that for $ = (0, 0)^ G the eigenvalues are 


Ai,2(<f>) = ^|[(l + l)-(2+(l^l)ni)0]Mi($) + (l±l)M2($)} 

_ fM2($) - 0Ml(<f>), 

1 (1 - (1 +ni)0)Mi($), 


and for $ = (0, 0)^ G the eigenvalues are 

Ai,2(<f>) = ^|(1 + 1 )Mi(<F) + [(1 ± 1) - (2 + (1 ± l)n2)0]M2(<f')| 

^ f(l - (1 +n2)0)M2($), 

- 0ti2(<f>). 


Within these eigenvalue characterizations we can distinguish the two types 


Aa(<F) := (1 - (1+ ni)0)M0$), for G 9^, z = 1, 2, 

Ab(<F) := %_*($) - 0zzi($), for <F G 9^, z = 1,2, 
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where the subindex 3 — i denominates the eomplementary index. The notation 
Aa(<h), A 6 (<I>) with new subindiees instead of Ai(<h), A 2 (<h) is used sinee the order 

Aa($) = Ai(<h) < A 2 ($) = Afe(<h) (35) 

eannot be always guaranteed. However, the notation Ai, A 2 is used whenever the 
order ( [35] ) ean be assured, e.g. at the origin, where 

Al(0) = Aa(0) = Uoo 2 < 'l^ool = Afe(O) = A 2 ( 0 ). 

Generally, order ehanges ean oeeur on both axes, as ean be seen from Fig.[^(b). 
The order ( [35] ) holds only in the absenee of sueh order ehanges. However, the oe- 
eurrenee of eoineidenee points does not affeet the subsequent elassifieation. In the 
sequel, it is shown how the eigenvalues Aa(<h) and Xb{^) assume the role of the 
eigenvalue of the first or the seeond family in dependenee of the value 0 . 

4.1. Inflection curves 

The highly non-linear strueture of the flux funetion makes it diffieult, if not 
impossible to obtain an explieit formula for infleetion manifolds liT^ . However, 
a eharaeterization for states on the eoordinate axes of the phase spaee ean be 
obtained, sinee the infleetion points there eorrespond to those of sealar equations. 
Namely, the eigenvalue Aa(<h) has the same form as the first derivative of the sealar 
flux funetion (N = 1), whereas there is no eorrespondenee for the additional 
eigenvalue Ab(<h), whieh emerges for the 2 x 2—system (N = 2). The derivative 

+ i = l,2 (36) 

is well defined for 0 e [0, 1). From ( |36| ) we have that A^0) may vanish at 0 = 1 
and at 0 ^ = 2/(1 + rii); the first ease is less relevant beeause it represents a state 
on a vertex of the phase spaee. We have that 0^ G (0, 1) sinee n* is larger than 
one due to eondition (S2). Thus the state (jf represents an infleetion point and the 
states <!)/ sueh that 0 j = 0 ^ and 0 ^ = 0 for i 7 ^ j G { 1 , 2 } belong to an infleetion 
manifold. The pertinenee of the infleetion points to the first eharaeteristie family 
ean be shown for general parameter settings. 

Lemma 7. There is at least one inflection point on each axis of the phase space 
0. having locations := (2/(1 + tt-i), O)^ and $2 •= (O5 2/(1 + ^2))^. 
Moreover, the corresponding eigenvalues belong to the first characteristic family. 
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Proof. The location of the inflection points on the axes is obtained by finding the 
zeros of the corresponding eigenvalue derivatives (j^. 

The attribution of the inflection points to a certain characteristic family is done 
by comparing the magnitudes of the eigenvalues Aa($*) and Afe(<f)^), both given 
by Substituting 0* = 2/(1 + rij) one observes that 


A„(<f>/) + = (1 - n,ctf)um) = 


l-rij 

l+Ui 




Since Ui, Uj are nonnegative and (1 — ni)/{l + nf) is always negative by the 
parameter setting (S2), we have that 

A,(«i>t)<«,(«i>/)-0V(«i>/) = A,($n. 

Thus, in the inflection points on the axes, the eigenvalue Xa corresponds to the 
first (smaller) eigenvalue Ai such that the eigenvalue order ([^ holds. □ 


4.2. Eigenvectors 

The Jacobian matrix ( [^ with components ( [3T| ) can be written as a rank two 
modification 


2 

J = D + Y^ akbl = D + BA^ (37) 

k=l 

of the diagonal matrix D = diag(ni, V 2 ) G where 

-0iA (1 m + 

\(j)2U2{^) -02/’ V + <1>'^W'(<1)) 

In a rank two modification one introduces column vectors 0^,6^ G ^ G 
{1,2}, in matrices A, B G which both have rank two. The formula for the 
components of a right eigenvector can be deduced from secular equations and is 
given as a function of eigenvalues ifTTl 



VijiX) 


1 Ll 





ajhl Y 

vk-x)y 


(39) 


where the parameters 

a] = 1, , 
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h] = h] = - fj, 


(40) 










correspond to the entries of the matriees A and -B in (381. Whereas formula (39) 


holds for general rank two modifieations of form ( [37] ), whieh are valid for system 
of arbitrary size, the following Lemma breaks it down for the eonsidered model 
with N = 2. 


Lemma 8. The right eigenvectors of the Jacobian matrix ( |30| ) are functions of the 
corresponding eigenvalues Ai(<h), i = 1, 2, given as 


^ U(ni(<h) - A,(«h)) - 


(41) 


Proof By the abbreviation 


A- = 

Vi — X 


one ean rewrite the eigenveetor ( [39l ) with parameters ( |40| ) eompaetly as 


1 + ^ 


k - Up 


k=l 


, J = l,2, 


r,,($,A($))=A, 

with the eomponents 

r*i = Ai [l + A 2 (m 2 - u'f)] = 
r..=A,[l + A,K-„a] = ^^_^ , 

Multiplying by (ui — A)(u 2 — A) the eigenveetors get the form ( |4T] ). 


+ 


0102(^1 - Uf 


V2-X (fi - A)(u2 - A) ’ 

0102(^2 - Ui) 


+ 


□ 


Speeial eigenveetors, in partieular those for the values of the Hugoniot loeus 
of the origin, ean be obtained from either exploiting the strueture of the Jaeobian 
matrix ( [SO] ) or from using the analytieal form of the eigenveetors ( |4T] ). For in- 
stanee, the eigenveetors on a set (3($*) that eorrespond to the seeond eigenvalue 
A 2 (‘h*) = ni($*) have the form 


r2($*) = 


(42) 


In the origin one has ri{0) = (0, I)"’" and r2(0) = (1, 0)"’". 
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4.3. Illustration of a benchmark example 

For illustration, a benchmark example (Example 1) is considered with param¬ 
eter setting 


r>ooi = 1, Voo 2 = 1/2, and ni = 4, ria = 3, (43) 

that satisfies specifications (S1)-(S3). The integral curves in the coordinate 
plane are shown in Fig. (left). The arrows point into the direction of increasing 
eigenvalues. The first characteristic family is crossing all lines of constant 0 and 
the second family is connecting the axes. It can be recognized that the second 
family is genuinely nonlinear, whereas genuine nonlinearity of the first family 
is lost at an inflection manifold, where the corresponding eigenvalues take their 
minimum. 

The direction of increasing eigenvalues of the second characteristic family 
switches at the contact manifold. This direction switch impacts on the solution 
structure of the Riemann problem RP(0, <F’''), depending on which side of the 
contact manifold the right state <!>■*■ is positioned, see Fig. (right) for results of 
simulations by a finite difference method. Only a slight change of the Riemann 
data provokes a fundamental change of the solution path in the phase space. 



Figure 1; Phase space of Example 1 with parameters Left; Integral curves in the 

0102—coordinate plane; First family (crossing all lines with constant 0) and second family (con¬ 
necting the axis). The integral curves point into the direction of increasing eigenvalues; contact 
manifold (dashed), inflection manifold of the first family (solid). Right; Solution of Riemann 
problems RP{0, (0.2, 0.25)), RP(0, (0.2,0.35)) by a finite difference method (low resolution). 
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5. Quasi-umbilic point 


The best examined cases of the loss of hyperbolicity are related to an “umbilic 
point” [[20ll^l29l . which is an isolated point where strict hyperbolicity fails. A 
description of new types of isolated points with loss of hyperbolicity is given in 
tn\ . In particular, if there is a connected set of points with loss of hyperbolicity, 
then more refined characterizations are needed. For instance, a generalization of 
the umbilic point is the coincidence point. 

Definition 6. A point $ is called a coincidence point of the PDE Q with flux 
function /($) if the eigenvalues of the Jacobian matrix J{^) of the flux function 
coincide at this point. We say that a coincidence point is an umbilic point if it 
satisfies the following conditions: 

(HI) The Jacobian matrix J(^*) is diagonalizable. 

(H2) There is a neighborhood V of^* such that J(^) has distinct eigenvalues 
for all 

Typically, we would expect that umbilic points are isolated coincidence points. 
However, the previous definition of an umbilic point may fail in either of the 
two conditions. Following [|^ . we classify an isolated coincidence point where 
condition (HI) fails as a quasi-umbilic point. It seems that the classification of 
quasi-umbilic points appears at the first time in [l26l . in IITtII a detailed description 
is given. 

Lemma 9. There is a unique coincidence point on $ G d^. Upon the choice of 
parameters there may exist up to two coincidence points on $ G d"^. 

Proof By definition, a state $ is a coincidence point if and only if Aa($) = Ab(<I>) 
holds. Namely, from the eigenvalue expression ( [3^ , one has a coincident point 
on axis $ G if there exists some fi G [0, 1] such that (1 — riiffuifb) = 
holds. Therefore, we must look for such a value fi which is a zero root of 

Ri{f) - V 002 /V 001 = 0, where 

Notice that by (SI) and (S2) the ratio V002/V001 is smaller than one and that the 
power oi—77-2 is positive. Since i?i( 0 ) = l,i?i(l/ni) = 0 and (</>) = —ni(l — 
0)"'i-"'2 _ _ 77,2)(1 — nicffjfi — is negative for all f G [ 0 , 1 /ni], 

we conclude that Ri{(j)) = V 002 /V 001 occurs at a single point, say 0“, which is 
the unique zero root of i?i( 0 ) — V002/V001 = 0 . Let us denote such a state as 
Qi ■■= iff Of. 
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Similarly, for the axis $ G 5^, from the eigenvalue expression ( [34| ), we look 
for a value 02 G [0, 1] (while 0i = 0) whieh is a zero root of R 2 { 4 >) ~ Vooi/voo 2 = 
0, where 

/?2(0):=(l-n20)(l-0)"S 0^1. 

Notiee that -R2(0) = 1 and -R 2 (l/'^ 2 ) = 0 hold, and sinee ni > n 2 > 1, in the 
limit lim 0 _>i_ R2{4>) = —oo holds. Moreover, i?2(0) is positive for 0 e [0, l/n 2 ) 
and negative for 0 G (l/n 2 , 1). As the ratio Vaoi/vao 2 is larger than one, there is 
no eoineidenee point if i?2(0) is always smaller than sueh a ratio. For 0 G [0, 1) 
we have 


i?'(0) = _n2(l-0)"^-"i-(ni-n2)(l-n20)(l-0)"^-' 
= [- 77 , 2(1 - 0) + {n 2 - 7ri)(l - 77,20)](1 - 


Thus, a single extremum of i?(0) oeeurs at 

^1 - 2^2 


( 77 l - -712 - 1)^2 


(44) 


There is one single eoineidenee point on the axis if i?2(0™) = 'i^ooi/'i^oo 2 - The 
existenee of two eoineidenee points oeeurs if i? 2 ( 0 ™) is larger than Vooi/voo 2 ', say 
Q2 := (0, Q3 := (0, 0^)"^ where i?2(02) = ^ 2 ( 03 ) = ^^ooi/^^oo 2 - □ 


The neeessary and suffieient eonditions for eoineidenee points on 0^ are dif- 
fieult to evaluate. Let us define the ratio W := Vaoi/voo 2 - On the one hand we 
know that V exeeeds one due to eondition (SI). On the other hand the powers 
77 ^ 1 , n 2 are also free parameters. Some ealeulations are possible if the differenee 
between ni and 772 is a natural number. For example, for 77i — 772 = 1, there are 
two eoineidenee points if 7r2 > W; for ni — n 2 = 2, there are two eoineidenee 
points if 77,2 > 2{W + a/IF^ — W); for rii — n 2 = 3, there are two eoineidenee 
points if 477,| — 271F (712 — 1)^ > 0. The inequality i?2(0m) < 1 holds for powers 
satisfying 772 + 1 < tzi < 2772 , and for 0m definied in ( |44l ). In this ease there are 
not eoineidenees on the axis 0^. 

For F = (0, 0)^ the Jaeobian matrix is an upper triangular matrix and for 
$ = (0, 0)^ a lower triangular matrix. In both oases the Jaeobian matrix is di- 
agonalizable. In view of Definition (|^, eondition (H2) is satisfied, but (HI) not 
neoessarily. Given that the diagonalization is valid for all <l> G 0^, the eoineidenee 
point Qi = (0“, 0)^ is a quasi-umbilio point, where the Jaeobian matrix has the 
form 
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is a Jordan block with * 7 ^ 0 and therefore not diagonalizable. This means that 
(HI) is violated, but (H2) is satisfied. Therefore, from Lemma we have the 
following result. 

Corollary 1. Any isolated coincidence point on the boundary turns out to be 
quasi-umbilic. 


(a) Bidisperse, example 1 



(b) Bidisperse, example 2 (c) Three-phase flow ELD 




Figure 2: Extended domains. Concentration and saturation triangles are marked by a bold solid 
line, shaded regions contain states with complex eigenvalues of the flux Jacobian. In (a) the 
Example 1 reveals two quasi-umbilic points at the boundary of elliptic regions; the point Qi was 
detected by our analysis. In (b) the Example 2 shows the three quasi-umbilic points detected by 
the previous analysis, another one appears at the boundary of the elliptic regions. The plot in (c) 
represents the ELD where the two quasi-umbilic points also belong to the boundary of elliptic 
regions. 


For illustration of the behavior of the eoineidenee points the flux funetions @ 
are eontinuously extended beyond the phase spaee by eonsidering relative veloei- 
ties @ without the eut-off (|7]). The visualization in Fig.|^shows that the identified 
quasi-umbilie point happens to belong to an elliptie/hyperbolie boundary; thus, it 
is elear that Ai(Qi) = \ 2 {Qi) holds. Sueh a eoneept is introdueed in [[2^ . and 
here, their flux functions are also extended by letting the dominated quadratic 
terms to act outside the physieal domain. In this ease the quasi-umbilie points 
also belong to an elliptie/hyperbolie boundary, see Fig.[^(e). 

For the bidisperse model with parameter specifications (S1)-(S3) two general 
examples are considered. Besides Example 1 with parameter setting ( |43| ), Exam¬ 
ple 2 has the parameters 


vioo = 1,^200 = 1/2, ni = 4.6, ^2 = 1.5. (45) 

For the ELD model (Equal-Lighter-Density fluids) ease in Rodrlguez-Bermudez 
and Marehesin (2013) the parameters are ehosen as pi = 2.0 > 1.0 = p2 = Pa¬ 
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6. Shock classification of Hugoniot locus of origin 

In this Section the types of shocks which are connected to the origin are classi¬ 
fied. The shock classification of Definition [^distinguishes three different types of 
admissible shocks, namely 1-Lax, 2-Lax and overcompressive shocks. The loca¬ 
tions of the shocks on the Hugoniot 'H{0) of the origin O = (0, 0)^ are identified 
in LemmaQ it includes in particular two contact manifolds. The shock classifica¬ 
tion builds on the eigenvalue analysis of Section j^by comparing the shock speed 
and the eigenvalues in each state on the Hugoniot locus. 

A key feature for the determination of the shock type consists of the order 
switch of the relative velocities at the threshold concentration cfy^, 

rMi($) > for 0 e[O,0 '"), 

< Mi($) = M 2 ($) for (46) 

< M 2 ($) for 

which is a direct consequence of the definition of the relative velocity 0, o, 
together with the specifications (SI) and (S2). 

A visual guide for the shock classification of shocks between the origin and 
states on the edges and 0^ is given by Fig. [sj where the shock speeds are com¬ 
pared to the characteristic speeds. The shock classification is essentially obtained 
by speed comparisons. In the following Theorem, the shock characterization is 
established for general parameter choices. 

Theorem 2. The Riemann problem 

RP(0, $), $ e H{0) = d^Ud‘^U C($") U C(<I)'^) 

connecting the origin to a state of its Hugoniot locus is solved by a single shock 
of speed a = (j(0, $), which can be classified as follows according to Definition 

m 

Classification for ^ = ((^ijO)"*^ G d^: 


2-Lax: 

Ai,2(<^>) < 

a < \ 2 {P) 

and Ai(0) < a 

for 

01 

e 

(0, c 

r). 


Ai,2(‘^*) < 

a < \ 2 {P) 

and Ai(0) = a 

for 

01 

= 



OC: 

Ai,2 

($) < a < 

Ai,2(0) 

for 

01 

e 

(0", 

0"). 

1 -Lax: 

Ai(<l>) < cr 

< Ai,2(0) 

and a = A 2 ($) 

for 

01 

e 

{0", 



Ai(<l>) < a 

< Ai,2(0) 

and a < A 2 (‘l’) 

for 

01 

e 

(0^ 

0“). 


(47) 
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Figure 3; Eigenvalues and shock speeds of shocks from the origin to (a) the states $ = (0, (f)) G 
9^, (j) G [0, 1], (b) the states <i) = {(j), 0) G 9^, (j) G [0, 1], Bold solid curves represent the 
shock speed cr(0, $). The eigenvalues (characteristic speeds) given in ( |34l i are represented as 
light solid curves for Aa($) and dashed curves for A;,(<1>). The dotted horizontal lines represent 
heights corresponding to cj)^, and </>“, ordered from top to bottom. 


The threshold value 0 °^ is given as 

0“^ = 1 - V^oo2/'yool- (48) 

Classification for <l> = (0, 4 ) 2 )^ £ 


1 -Lax: Ai($) < a < Ai,2(0) and a < A2($) 
Ai($) < a < Ai,2(0) and a = A2($) 
OC: Ai, 2 (<^>) < cr < Ai,2(0) 

Classification for ^ G 


for 02 e (0,0^), 
for 02 e {0^,0°°}. 
for 02 e (0^, 4 °°). 


1 -Lax: Ai($) < a < Ai_2(0) and a = A 2 ($). (50) 

Classification for ^ G C($°°).' 

1 -Lax: Ai($) = a < Ai^2(0) and a = A 2 ($). (51) 


Strict inequalities in the shock type classification ([77])-([5T]) mean that the cor¬ 
responding shocks are not characteristic shocks. For instance, for <1) G 0^ there 
are no 1-Lax shocks characteristic at the left datum O. 
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Proof. The proof is done in two main steps: 

1. - The first step is to relate the eharaeteristie speeds Ai,2(0) and Aa,fe($) to the 
shoek speed cr(0, $) in dependenee of the position of 

2. - The seeond step eonsists in determining whieh shoek elassifieation applies ae- 
eording to Definition 

Table gives an overview on the speed magnitudes for right states on the axes 
and d"^. 



< cr(0, $) 

> cr(0,<f)) 


8 ^ : 0 < 0°" 

A„,„(<f)) Ai(0) 

A2(0) 

2-Lax 

(r,r) 

Aa,6(<^>) 

Ai,2(0) 

OC 

0 > 0^ 

Aa(<^>) 

Afe($) Ai,2(0) 

1-Lax 

8 ^ \ (f) < 

Aa(<^>) 

Afe($) Ai,2(0) 

1-Lax 

0 > 0^ 

Aa,6(<^>) 

Ai,2(0) 

OC 


Table 1; Shock classification on Hugoniot locus of origin for states on the axes and d^. First 
column; Location on the axis. Second column; Characteristic speeds less than cr(0, $). Third 
column; Characteristic speeds larger than cr(0, $). Forth column; Shock type. 


A 2 ( 0 ) eompared to cr((9, $) on and 

On the axes G d\0 the shoek speed is limited as 


cr(0,<f>) = ni($) < Vooi 


Vool = A2(0), 

Voo2 = Ai(0) < A2(0), 


if i = 1, 
if i = 2, 


(52) 


sinee the shoek speed 


cr(0, $) = = (1 - </>*)«*($) = ^00^(1 - fiTf * = 1, 2, 

is monotonieally deereasing on both axes and 8“^. (Note that ( |52l ) exeludes 2- 
Lax shoeks on edge 8“^.) 


Aa($) eompared to cr(0, $) on 8^ and 8"^ 


24 











The eigenvalue Aa(<l)) as speeified in ( [34| ) leads to eommon properties for both 
axes. For all points G i G {1,2}, along the axes one has 


Aa(<l>) = [!-(! + ni)(f)i]ui{^) 

< (1 - (/>*)«*($) = ni($) = cr(0, $), i = 1, 2, 


(53) 


due to properties that apply on the axes, namely, the speed eonvertibility ( [291 ) and 
the shoek speed ( |27| ). 


A6($) eompared to a{0, $) on and 

The eomposition of the eigenvalue Afc($) depends on the inequality ( j^ be¬ 
tween the relative veloeities Mi(<F) and M 2 ($), and induees quite different behav¬ 
iors on the edges. On the axis the eigenvalue A6(<f>) behaves as 


r A6($) < a(0, <F) for G (0, 0"), 

I A6($) = a(0, <F) for 0le{0^0~}, (54) 

[\($) >a(0, <F) for 0lG(0^0~), 


beeause on d^, the inequalities 
interval 0i G [0, cf)^) the eigenvalue (34) satisfies 


between mi(<I') and M 2 (* 1 >) implies that for the 


Afe($) = M 2 (<F) - 0iMi($) < (1 - 0i)mi($) = cr(0, <!)), 


and for the other interval, where 0i G (0^^, 0“], the inequality sign switehes. 

On the axis 0^ one has 

r Afe($) > a(0, <!)) for 02 G [0,0^), 

< A6($) = cr(0, <!)) for 02 G {0=^,0“}, (55) 

[a 6($) < cr(0, <!)) for 02G(0^,0“), 

sinee, from the inequalities ( |4^ for 02 G [0, cfy^), one obtains 

A6($) = Ui(<f)) - 02M2(<1>) > (1 - 02)^i2(*l’) = <^{ 0 , $) 


and for 0 G (0'', 0°°] the inequality sign switehes. 

On the edge 0^ the value of Afe($) determines whether the shoek is 1-Lax 
or over-eompressive. The shoek type is deeided by the relative magnitudes of 
Mi($) eontra M 2 (*fi): The equality Mi($) = M2(*1*) only holds for 02 = (f)^ (and 
02 = 0“)- On edge there is a eoineidenee of eigenvalues, see also Lemma|^ 


25 






Ai(0) compared to a{0, $) on and d'^ 

Note that, along both axes a(0, $) is a monotonically decreasing function 


with (t( 0, $“) = 0. By ( [5^ is assured that aiO, $) < Ai(0) holds for all states 
on the axis d‘^\0. However, this does not hold for all states on the axis d^\0. 

Since lim£_ 5 .o a(0, (e, 0)) = A2(0) > Ai(0) > 0 on the axis there exists a 
state = (0°^, 0)^ such that (t( 0, = Ai(0). In this state the equality 

a(0, = Uooi(l - (/>ar' = Voc2 = Ai(0) 


holds, from which the characterization of 0^^ by ( |4^ can be deduced. The shock 
speed ct(0, <h) relates to Ai(0) in dependence of 0°^ as 


f cr(0, <!)) > Ai(0) for 01 e [0, 0“^), 

I a(0, <h) = Ai(0) for 0i = 0", (56) 

[cr(0, $) < Ai(0) for 01 e (0*^, 0°°]. 


Together with the previously established ( p^ , ( |5^ and ( [54| ) the discrimination 
implies that, on the axis d^, for 0 < there is a 2-Lax shock, whereas for 
<1)'^ < 0 < 0"^ the shock is over-compressive. 

By ( [5^ , ( |54l ), for 0 G (0^, 0“), on axis there is a clear separation between 
the shock speeds: 


Ai(<h) = Aa(<I)) < a(0, $) < Ab($) = A2(<I>). 


(57) 


This separation allows an association of the eigenvalues according to ( |35| ). The 
same way, by ( |5^ , ( |55] ) for 0 G (0,0^) on there is a clear separation be¬ 
tween the shock speeds ( [57| ) that establishes the eigenvalues order ( [35] ). For 
0 G (0'^,0°°) there is no clear eigenvalue separation, which however does not 
affect the fact that the shock is over-compressive. 

Now we are able to conclude the shock classification on the axes: Putting the 
inequalities ([52l), ([5^, ([55]) together gives the shock classification (|49|) for right 


states on the axis d^. Putting the inequalities ( |52j ), ( |5^ , ( |54| ), ( |5^ together gives 
the shock classification ([47]) for right states on the axis 


Contact manifold (7(<h'") 

The shock classification for states on the contact manifold $ G (7(<I)’^) assures 
that all states are connected to the origin by a 1-Lax shock: Since the functions 
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and n 2 (*^*) take their maximum in the origin O and the origin is not part of 
the contact manifold, one has the strict inequality 


a{0, $) = < Uooi = Ai(0) < Voo2 = A2(0). 


(58) 


By the eigenvalue characterization ( [^ for states on and noting that A$ > 

0 for all <h G one gets 


Ai($) < a{0, <h) = ui($) = A2(<h) 


(59) 


for all <h G C(<h''). The properties ( |58] ) and ( |59l ) are summarized in the shock 
classification i 


Contact manifold d°° = C(<h“) 

For states on the maximum packing manifold $ G d°° one has vanishing 
eigenvalues, 

0 = Ai($) = A2($) = a{0, <h) < Ai(0) < A2(0), 

which leads to classification ( [5T] ). □ 





Figure 4; Different shock types for states connected to the origin. Continuous curves are 2-Lax 
shocks, dashed curves are 1-Lax shocks and dotted-dashed are over-compressive shocks. The 
Hugoniot locus of the origin comprises the three edges and the contact manifolds C(<i)’^), C(<i)°°). 
Also some 2-Lax shock curves from arbitrarily chosen generic points and are 

indicated. Reference values on the horizontal axis are (j)'^, (j)^ and (p°°, on the vertical axis 
are (j)^ and 0°°. 
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Figure [^displays the different shock types for states connected with the origin, 
as elaborated in Theorem]^ The 1-Lax states of 'H(O) are shown as dashed lines. 

States on the contact manifold satisfy the 1-Lax conditions being char¬ 
acteristic in the second family at the right state. In the same manner, the states 
on the maximum packing = d°° satisfy to be 1-Lax shocks which are 

characteristic at the right state with respect to both families. On the edge d'^, at 
<I)^ = (0, and = (0, 0°°)"'" the shocks are 1-Lax and characteristic at 
the right states, say cr(0, <I)^) = A 2 (<I’^) and a{0, = A 2 (<F“). On the edge 

d^, at 0)"'", <1)°' = (0'^, 0)"'" and = (0“, 0)"*" the shocks are 2- 

Lax, over-compressive and 1-Lax, respectively, being characteristic at left or right 
states, say a{0, <I>^) = Ai(0), a(0, = A 2 (<F'^) and ct(0, $“) = A 2 (<F“). 

Several states located in Fig.[^(b) have characteristic shocks: If $ G then 

the shock is right characteristic for the second family, if 0 = 0“ then it is char¬ 
acteristic for both families, and if $ = (0, 0^^)^ then the shock is 2-Lax but left 
characteristic in the first family. 

7. Riemann problems 

In this Section, we construct the solution of the Riemann problems RP(0, $) 
and RP($, where O = (0, 0)"'", G and <I) G "D^oo is a generic right 
or left state in the phase space. These Riemann problems are derived from the 
standard initial condition ( [T?] ). 

7.1. The Riemann problem HF {O, 

The behavior of solutions to the Riemann problem RP(0, <!)) depends on the 
position of <I) in the interior of with respect to the contact manifold 
which splits the domain into two regions: 

V~ := {$ G Vq,oc : 0 < 0^}, 

p+ ;= {$ G P$oc : 0 > 0"}, 

where C(<!)'') = cl('D~) fl cl(P0) marks the intersection of the closures of those 
two regions. The main difference between both domains is the opposite direction 
of the characteristic speeds on the integral curves. In V~ the second eigenvalue 
increases from the edge to the edge d^, and, reversely, in the second eigen¬ 
value increases from the edge 0^ to the edge d^. See the orientation of the second 
rarefactions in Fig. where the arrows point into the directions of increasing 
eigenvalues. 
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A Riemann solution from O to any state $ in generally eonsists of a 1- 
Lax shook followed by a 2-Lax shook. For belonging to the middle state that 
interseots the 1-wave with the 2-wave is denoted by , and for belonging to 
the middle state is denoted by . In both oases the middle state belongs 
to the 1-Lax loons of the origin O. Sinoe both waves are shooks any middle state 
or is looated at U 'H(O). It turns out that belongs to the edge 
and belongs to the edge d^. 

For a right state <I) G the solution of RP{0, $) eomprises the 1-Lax 
shook from O to a state <I>^ = (0, 0^)^ sueh that 0^ G (0, 0"^) and the 2- 
Lax shook from <I>^ to <I). Similarly, for a right state <I) G T>^ the solution of 
RP{0, $) eomprises the 1-Lax shook from O to a state = (0^, 0)^ sueh 
that 0^ G 0°°) and the 2-Lax shook from to <h. Sueh oases are depleted 
in Fig.Q the former as $_ G V~ and the latter as $+ G 

The solution of the Riemann problem RP{0, $) with <I) on the Hugoniot loous 
of O eomprises a single shook from O to $ with a olassifioation that depends on 
the position of $ and is elaborated in Theorem 

1. For <I> = ( 0 , 0)"'" with 0 G ( 0 ^^, (j)^) the shook is over-eompressive, 

2. For <I> = (0, 0 )"'" with 0 G ( 0 ^^, 0 °°) the shook is over-eompressive, 

3. For $ = ( 0 , 0)"'" with 0 G (0, 0 °^) the shook is 2-Lax, 

4. For any other <h in 'H(O) the shook is 1-Lax. 

For a state <h G C (<!)''), the solution of the Riemann problem RP{0, <h) eon¬ 
sists of a single 1-Lax shook from O to <I). The two oonseoutive shoeks from O 
to (0, 0^)^ and from (0, 0^)"’" to ^ have both the same speed ui(0^) = U 2 ( 0 ’^), 
whioh in turn is the same speed of the direot shook from O to <h. Sinoe all these 
shooks have the same speed, the middle state (0, 0"")^ is “invisible” in the solution 
profile in the physioal spaee; this is beoause of LemmaTherefore, the solution 
strueture for <I) G "D^oo is represented as 

o -)■ $ -> $. 

For <h G C(<I)^) the middle state <I>^ may assume any other state on the same 
oontaot manifold and even oollapse with <h. 

Please refer to Fig.|^and notiee that as any $_ G V~ tends to a $ G the 
middle state tends to (0, 0^^)^. Similarly, as <h+ G tends to a <I) G C(<I)’^), 
the state tends to (0^^, 0)^. Thus, we notiee a oontinuous dependenoe of 
solutions to the Riemann problem RP(0, <!)) on the right datum. 
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7.2. The Riemann problem RP{^, ^°°) 


In this Section, the dependence of the solution structure on the exponents ni 
and 77-2 is illustrated by two examples, Example 1 and Example 2 with the cor¬ 
responding parameter setting ( |4^ and d^ , respectively. In the Example 1, the 
solution consists of a simple 1-wave comprising a shock followed by a rarefaction 
or a single rarefaction. The solution structure of Example 2 depends on the exis¬ 
tence of one or two detached inflection curves for the first characteristic family. 
As is a contact manifold, the goal is to consider a generic state $ G 

from which waves reaching any state in G d°° are constructed. 


Construction for Example 1. The construction for Example 1 can be orientated 
by the integral curves and the inflection manifold in Eigure[^ In Eigure[2(a) the 
inflection manifold Xi of the first family is shown as a solid curve, which is almost 
a straight line, connecting the states = (0.4, 0)"'" and = (0, 0.5)"’", compare 
Eemma|^ Therefore, a wave curve starting at a state on the upper-right hand side 
of Ii follows a centered rarefaction of the first family until the maximum packing 
manifold is reached. The final rarefaction point is (0, 1)^, unless the starting 
state is on the axis, in such a case, the final rarefaction point is (1, 0)^. 

The characteristic velocities near the contact manifold C(<1)°°) are close to 
zero, and the characteristic directions are close to r 2 = ri = (1, —1)^, see ( |4^ . 
If the left state $ is on the lower-left side of the inflection Xi then the wave is 
obtained by a backward 1-wave construction. Namely, all shocks from a state 4) 
at the lower-left hand side of Xi are connected to a state 4)^ G on the right 
of Xi satisfying 4)^) = Ai(<l>^). Therefore, the solution for such a state 4) 
comprises the 1-Eax shock from 4) to 4)^ which is characteristic at <1>^, and the 
rarefaction curve from 4)^ to (0, 1)^, or (1, 0)^ if $ belongs to . 


Construction of Example 2. The construction of Example 2 is visualized in Eig- 
ure 1^ The integral curves of the first family change their growth direction at two 
detached parts of the inflection manifold Xi, namely T and B, see Eigure[^(a). 
This growth direction change is the key in the distinct structure of the solution of 
RPifb, solutions in both Examples. 

In Eig. (b) the light and dark shaded regions represent right states $ for 
which the solutions are similar to that of Example 1: 

1. When the states belong to the upper comer above T inside the light shaded 
region, then the wave curve comprises a single 1-rarefaction connecting the 
state 4) moving along increasing eigenvalues Ai towards = (1, 0)"""; 
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(a) 1-Rarefactions (b) Structure in the triangle (c) Typical solutions 



Figure 5: Curves of Example 2. In (a) the first family rarefaction curves are plotted as continuous 
curves, the arrows point into the direction of increasing eigenvalues. The dash-dotted curves are 
the inflection manifold X\ splitted into the two branches T and B. In (b) the solid curve Ci i is 
the double-contact of the first family (notice the two components, one in the light region and one 
another in the white region); A is the state where a limit shock curve from T does not cross B. In 
the presentation of the wave curve solutions in (c) the continuous curves correspond to rarefaction 
fans and the dashed curves correspond to shock waves. 


2. For a left state $ in the dark shaded region, a eharaeteristie 1-Lax shook to 
a middle state in the light shaded region erosses T, from this middle state a 
first family rarefaetion follows to the maximum paekage eoneentration. 

The eonstruotion of $“) solutions for in the white region of Fig¬ 

ure 1^ (b) is outlined in the sequel. There are three eases of right eharaeteristie 
1-Lax shoeks that oonneet a state $ to a state G 'H{^) suoh that = 

Ai(<h^) : 

(1) The shook ourve erosses onoe 

(2) The shook ourve erosses onoe T ; 

(3) The shook ourve erosses twioe B and onoe T. 

In Figure[^(o) the left states of the oases (1) are represented by $1 and $ 2 , the 
ease (2) is represented by point 2 above and, ease (3) is represented by < 1 ) 3 . The 
wave ourves starting from $ 1 , ^2 and $3 have the following strueture: 

I 1—Lax liM 1—rar iR| 1—Lax I i (7 1—rar ino 

$1 -). -). -). 

I 1—Lax I I 1—rar iR| 1—Lax I i (7 1—rar ^no 

^2 -^ \^2 -^ 1 -^ 1*^2 -^ ) 

I 1—Lax liM 1—rar ino 

$3 -). -i. 

where the symbol I indieates where a shook is eharaeteristie. For ease (3) the 


(60) 

(61) 

(62) 
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construction of the shock curves proeeeds as before, namely a 1-Lax shoek fol¬ 
lowed by a 1 -rarefaetion eonneeting a middle state in the light shade region 
in Fig. (b) with For ease (i), with z = 1 or 2, the eharaeteristie shoek 
($j, preeedes a 1 -rarefaetion from towards S at a state <Ff on Ci^i, from 
there another 1-Lax left and right eharaeteristie shoek eonneets to a state <f)f on 
the other side of Ci i. From the wave eurve is terminated by the 1-rarefaetion 
to 

If the state < 1)2 is approximated to a state < 1)3 then a eontinuous ehange in the 
solution profile is observed sueh that the wave group ( [M| ) eollapses to the wave 
group (j 6 ^. Indeed, if $2 comes eloser to < 1 ) 3 , then eomes eloser to <I>f and 
^2 to ^ 3 ^ ill sueh a way that the shoek speeds a{^ 2 , and cT(<I)f, ^ 2 ) ap¬ 
proximate cr(<I) 3 , <I)f^), until, in the limit, the 1 -rarefaetion wave from to <I)f 
disappears. 
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